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ABSTRACT 

Time-dependent insolation in a planetary atmosphere induces a mass quadrupole upon which the 
stellar tidal acceleration can exert a force. This "thermal tide" force can give rise to secular torques 
on the planet and orbit as well as radial forces causing eccentricity evolution. We apply this idea to 
the close-in gas giant exoplanets ( "hot Jupiters" ) . The response of radiative atmospheres is computed 
in a hydrostatic model which treats the insolation as a time-dependent heat source, and solves for 
thermal radiation using flux-limited diffusion. Fully nonlinear numerical simulations are compared 
to solutions of the linearized equations, as well as analytic approximations, all of which are in good 
agreement. We find generically that thermal tide density perturbations lead the semi-diurnal forcing. 
As a result thermal tides can generate asynchronous spin and eccentricity. Applying our calculations 
to the hot Jupiters, we find the following results: (1) Departure from synchronous spin is significant 
for hot Jupiters, and increases with orbital period. (2) Ongoing gravitational tidal dissipation in spin 
equilibrium leads to steady-state internal heating rates up to ~ 10 28 erg s _1 . If deposited sufficiently 
deep, these heating rates may explain the anomalously large radii of many hot Jupiters in terms of 
a "tidal main sequence" where cooling balances tidal heating. At fixed stellar type, planet mass and 
tidal Q, planetary radius increases strongly toward the star inside orbital periods < 2 weeks. (3) There 
exists a narrow window in orbital period where small eccentricities, e, grow exponentially with a large 
rate. This window may explain the ~ 1/4 of hot Jupiters which should have been circularized by the 
gravitational tide long ago, but are observed to have significant nonzero e. Conversely, outside this 
window, the thermal and gravitational tide both act to damp e, complicating the ability to constrain 
the planet's tidal Q. 

Subject headings: planets - tides 



1. INTRODUCTION 

A number of puzzles have arisen for the gas giant ex- 
C*") . oplanets orbiting close to their parent stars, the "hot 
C***~ ' Jupiters." 

A large fraction of transiting hot Jupiters are observed 
y—i • to have radii far larger than the radius of Jupiter, Rj, 
' implyi ng high temperatur es deep in the planetary inte- 
Q\ \ rior (|Burrows et al.ll200Cf ). Strong irradiation has been 
■ found to have an insulating effect on planets, s lowing 
^ | their cooling and contraction (jBurrows et al.ll2000| ). This 
• i-H , effect can explain radii R p ~ (1.0 — 1.2) x Rj, but is in- 
sufficient to explain the radii of a significant fraction of 
the p opulation with R p ~ (1.2 — 1.8) x Rj (see iFortnevI 
?3 i |2008l for a recent review). A powerful internal heat source 
must be acting to prevent these planets from contracting. 

Circularization of the planet's orbit by dissipation of 
the gravitational tide has been invoked to explain the 
small eccentricities of m ost hot Jupiters (Rasi o et al.l 
119961: iMarcv et al1ll997ft. Using Ju piter's inferred tidal 
Q = 10 5 - 10 6 (|Yoder fc Pealdll981[ L orbits are expected 
to be circular out to orbital periods P OI b < 1 week. Curi- 
ously, a large fraction of the populati on that should have 
been circularized in ~Gyrs i s not (Do bbs-Dixon etlrtl 
12004 iMatsumura et al.l 120081 ). The zero and nonzero 
eccentricity planets occupy the same orbital period 
range. Suggested explanations are tidal interaction with 
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a rapidly rotating young star (|Dobbs-Dixon et ai¥ 2004) , 
perturbations fr om other planets, or a la rge range of tidal 
Q = 10 5 - 10 9 (|Matsumura et alj 12008]) . In the contin- 
ued absence of a detectable perturber with appropriate 
mass and orbit to account for the observed eccentric- 
ity, a mechanism involving an isolated star and planet is 
needed. It is also not clear why planets with presumably 
similar internal structure and orbits should have such 
different levels of internal friction. 

In the absence of direct observations, it is commonly 
assumed that the rotation of present day hot Jupiters is 
highly synchronous. This assumption is motivated by the 
short (~ Myr) synchronization time using Q = 10 5 — 10 6 
for the gas giant. However, despite the short expected 
synchronization time, hot Jupiters are particularly sus- 
ceptible to develop asynchronous rotation. Gas giants 
are relatively frictionless compared to terrestrial plan- 
ets with Q ~ 10 — 100. Even weak opposing external 
torques may compete with the gravitational tide to pro- 
duce asynchronous spin. 

The proximity of hot Jupiters to their parent star, 
~ 10 2 times closer than Solar System giants, implies 
an insolation stronger by ~ 10 4 and stellar tidal forces 
stronger by ~ 10 6 . Insolation and tidal forces may then 
play a far more important role for close-in planets. In 
this paper we discuss "thermal tide torques", which are 
created through the interplay of stellar irradiation and 
g ravitational tidal forc es. 

iGold &: Soterl (|1969h first applied thermal tide torques 
to explain the observed rotation rate of Venus. A num- 
ber of more detailed studies followed (e.g. Ingersoll and 
Dobrovolskis 1981). Venus resides in state of slow ret- 
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rograde rotation. Gravitational tides acting alone would 
have synchronized Venus' spin to high precision. Further- 
more, Venus' spin seems to be in resonance with Earth's 
orbit. Tides from Earth acting on permanent quadrupole 
in Venus should have a negligible effect in comparison to 
solar gravitational torques. The possibility that torques 
from solar gravitational and thermal tides nearly balance 
that allows a small effect such as tidal forces from Earth 
to have an observable effect. 

Previous studies of the thermal tide have assumed a 
thin atmosphere on top of a solid surface, as appropri- 
ate for Venus and Earth. In addition, they have as- 
sumed that the gravitational tide is strongly dissipative 
(Q ~ 10 — 100). In applying the thermal tide to close-in 
gas giant planets, a difference between our model and 
those for Venus and Earth is that we assume an opti- 
cally thick atmosphere, rather than ground, completely 
absorbs the stellar radiation. We include flux-limited ra- 
diative diffusion of heat above and below this absorbing 
layer. We find thermal diffusion is crucial in order to 
accurately estimate the magnitude of the thermal tide, 
particularly for long forcing periods. In applying our re- 
sults, we find that the thermal tide can affect not only 
the planet's spin, but also the eccentricity and radius for 
close-in planets. 

The plan of the paper is as follows. $2] contains a qual- 
itative discussion of why thermal tide effects are impor- 
tant. There, we argue that spin, eccentricity and global 
energetics of hot Jupiters are determined by the compe- 
tition between the thermal and gravitational tides. $3] 
contains a detailed treatment of atmospheric response to 
time-dependent insolation, and the resultant quadrupole 
moments. Readers wishing to ignore these technical de- 
tails should skip this section. Basic equations and ge- 
ometry are introduced in Q3.ll the heating function is 
expanded in a Fourier series in §3.2( temperature pro- 
files from the nonlinear simulations are discussed in £13.31 
quadrupole moments are computed in £)3.41 linear per- 
turbation equations are derived in ^3.5i and analytic so- 
lutions in the high and low frequency limits are presented 
in £ 13.61 and 13. 7\ respectively. Applications to hot Jupiter 
spin, radius and eccentricity are in the following sections. 
Equilibrium spin frequencies are computed in £jH <|5] con- 
tains tidal heating rates, and a discussion of the energet- 
ics required to power the observed radii. ^Sl contains 
eccentricity growth/decay rates due to both the thermal 
and gravitational tides. Our most detailed results for the 
simultaneous equilibria of spin rate, planetary radius and 
eccentricity are contained in £JTJ In <J5] we qualitatively 
discuss how thermal and gravitational tidal forces may 
drive differential rotation which, in itself, may lead to 
dissipation. Summary and conclusions are presented in 
fj9] The Appendices contains a derivation of the torque, 
orbital evolution rates, and tidal heating rates due to the 
presence of mass quadrupoles of arbitrary form, as well 
as the mass-radius relation and core luminosity. 

2. BASIC IDEA AND INITIAL ESTIMATES 

Tidal forces become increasingly important with de- 
creasing orbital radius. Acting alone, dissipation of the 
gravitational tide leads to synchronous spin and circular 
orbit for the planet. Thermal tide torques complicate 
this simple model, leading to qualitatively different evo- 
lution. A conceptual introduction and order of magni- 
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Fig. 1. — Geometry of the semi-diurnal (m = 2) thermal and 
gravitational tidal density perturbations. The planet rotates coun- 
terclockwise. In order for the torque to push the planet away from 
synchronous spin, the density perturbation must lead the line join- 
ing the planet and star. 

tude estimates are given below. 

2.1. Induced Quadrupoles 

First we summarize the gravitational tide. The stellar 
tidal force raises a bulge on the planet, quantified to 
lowest order by the quadrupole moment 

Here R p and M p are the radius and mass of the planet, 
M* is the mass of the star, a is the semi-major axis of 
the orbit, and h t ~ R p (M+/M p )(R p /a) 3 is the height of 
the tide. In the absence of friction within the planet, the 
gravitational tide does not cause secular evolution of the 
orbital elements or spin. The presence of friction causes 
high tide to lag maximum tidal acceleration at noon and 
midnight. For the case of asynchronous spin, this implies 
that the tidal bulge is misaligned with the line joining the 
planet and star, as pictured in figure [TJ The stellar tidal 
force then torques the planet by "pulling" on the projec- 
tion of the quadrupole moment Qvsrav) perpendicular to 
the star-planet line. This small lag in phase or time is 
conveniently parametrized by the tidal Q parameter (e.g. 
Goldreich and Soter 1965). The effective quadrupole mo- 
ment contributing to secular evolution is 




where a is the tidal forcing frequency and n is the orbital 
frequency. The factor a/n ensures the misalignment of 
the tidal bulge vanishes for synchronous and circular or- 
bits (see £|3.4p . The torque, N, on the planet due to this 

quadrupole moment scales as N ~ n 2 <2ifJ av ' ) . 

Jupiter's Q is thought to be in the range Q = 10 5 — 
10 6 based on ex pansion of Io's orbit to form the Laplace 
resonance (e.g. lYoder fe Peali |1981[ ). Several studies 
have also found that hot Jupiters should have Q = 10 5 — 
10 6 to explain t he circularization of t heir orbit s on Gyr 
timescales (e.g. lOgilvie fe Lml 12004 IWull2005f ). These 
studies assume that the gravitational tide is acting alone. 
Here, we show that constraints on Q of extrasolar planets 
are complicated by the influence of the thermal tide. 

Another possibility for tidal phenomena comes from in- 
tense time-dependent insolation. Time-dependence may 



3 



arise either through asynchronous rotation or orbital ec- 
centricity. For circular orbits and asynchronous spin, the 
diurnal forcing frequency is a = n — f2, where fi is the 
planet's spin frequency. There is also power at the higher 
harmonics a = m(n — f2), where to is the azimuthal 
wavenumber. Since the diurnal component does not gen- 
erate a quadrupole moment which can couple to the stel- 
lar tidal field, the semi-diurnal component (to = 2) dom- 
inates the torque on the planet. The magnitude of the 
thermal tide quadrupole is roughly (see M3. 61 and 13. 7|) 



~ (<rt th y 



(3) 



where AM is the mass of the atmospheric layer heated 
by time-dependent insolation, 



tth 



AMCpT 
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the exponent S ~ 0.5 — 1.0 and F+ is the flux at the sub- 
stellar point. The thermal time, tth-, is the time required 
for the layer of mass AM to absorb/emit its thermal 
energy. 

One of our principal goals is to show that for the at- 
mospheres of interest, the thermal tide bulge leads the 
sub-stellar point, opposite the behavior of the diss ipative 
gravitational tide. Ilngersoll fc Dobrovolsk is (1978) made 
similar arguments for Venus based on two assumptions. 
First, the temperature of the absorbing layer lags the 
time of maximum heating, similar to our daily experience 
on Earth. Second, that fluid elements remain at roughly 
constant pressure. Regions of lower temperature then 
imply higher density and vice versa, and the high density 
regions ("tidal bulges") lead the sub-stellar point. In our 
analysis, we adopt this simplifying "constant pressure ap- 
proximation." The horizontal readjustment of mass nec- 
essarily requires zonal winds. In this paper we do not 
solve f or the wind structure; see lDobrovolskis fc Ingersolll 
(1980) for the wind structure resulting from a similar cal- 
culation applied to Venus. 

Both the gravitational and thermal tidal bulges can be 
torqued by the stellar tidal field. The angular momentum 
evolves due to the net torque. In spin equilibrium, the 
thermal and gravitational tide torques balance, in general 
leading to an asynchronous spin state. We estimate the 
ratio of quadrupole moments to be (eq|2] and [3|) 
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This order of magnitude estimate shows that for param- 
eters characteristic of the hot Jupiters the thermal and 
gravitational tide effects can be comparable. Setting eq|5] 
equal to unity allows one to estimate the equilibrium spin 
frequency since a cx n — Q. 

2.2. Orbit, spin, and thermal equilibrium 

Gravitational tidal dissipation converts rotational and 
orbital energy into heat. When the gravitational tide 



alone acts, the end state is synchronous spin and a 
circular orbit. As torque equilibrium is realized in ~ 
Myr, over Gyr timescales the only source of tidal heat- 
ing is finite eccentricity. In the absence of an exter- 
nal perturbation continually pumping the eccentricity, 
tidal heating is limited by the initial energy reservoir in 
the orbit ~ GM+M p e 2 /2a, which may be comparable 
to planet's binding energ y for initially eccentric orbits 
(|Bodenheimer et all 12001ft . 

Planets with anomalously large radii but small ec- 
centricity, such as HD 209458b with e ~ 0.01 
(jLaughlin et al.l 12005ft . are difficult to explain with the 
gravitational tide alone since the current eccentricity im- 
plies negligible tidal heating. In order to explain the 
large radius with transient tidal heating one has to in- 
voke that heating ended recently, and we are observing 
the planet before it has had time to contract. Since con- 
traction times for radii (1.2 — 1.8) x Rj are much shorter 
than a Gyr, this argument requires that we observe these 
systems during a special time in their evolution. 

Inclusion of the thermal tide qualitatively changes this 
picture. A steady state may be reached in which thermal 
and gravitational tide effects balance one another. In this 
scenario, gravitational tidal dissipation continues to act 
in the equilibrium state, due to asynchronous spin and/or 
eccentricity. The thermal tide forces are not inherently 
dissipative. However differential rotation set up in the 
atmosphere by opposing torques at different depths may 
cause significant heating. 

Tangential forces applied to the thermal tide 
quadrupole torque the planet away from synchronous 
spin. In equilibrium, the thermal and gravitational tide 
torques on the planet are in balance, setting the equi- 
librium spin rate. The timescale to attain equilibrium 
spin is shorter than the Gyrs age for orbital periods up 
to several months. Radial forces from the thermal tide 
quadrupole act to alter the orbital eccentricity. If both 
the thermal and gravitational tides act to circularize the 
orbit, the equilibrium state is again zero eccentricity. 
However, the thermal tide can also pump eccentricity, 
opposing the gravitational tide. In that CM SC. cl balance 
of thermal and gravitational tides implies nonzero values 
for the equilibrium eccentricity. The timescale to attain 
equilibrium eccentricity is shorter than the age for plan- 
ets with orbital periods shorter than ~ 1 week. Lastly, in 
the absence of tidal heating, planets cool and contract. 
If tides deposit heat deep in the convective interior of the 
planet, a thermal equilibrium is possible in which heat- 
ing of the core is balanced by outward heat loss at the 
radiative-convective boundary. We find thermal equilib- 
rium is achieved on a timescale shorter than the age for 
orbital periods < 1 week (see section [7]). 

The thermal tide torque is exerted on a layer near the 
photosphere of the stellar radiation. The depth depen- 
dence of the gravitational tide torque is uncertain as the 
nature of gravitational tidal dissipation is complicated 
and not yet fully understood. For instance, turbulent 
viscosity damping would mainl y occur just below the 
radiative-convec tive boundary (|Goldreich fc Nicholson! 
Il977t IWul 12005ft , shear layers due t o "wave attractors" 
reside deep in the convecti ve core (|Qgilvie fc Linl[200l 
IGoodman fc Lacknerl [2008), and upward propagating 
waves genera ted at the tropopaus e may break above the 
photosphere (Ogil vie fc Linll2004ft . If the thermal and 
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gravitational torques are not exerted at the same depth, 
differential rotation is induced, distinct from the usual 
thermally-driven circulation patterns. For simplicity, we 
defer effects due to differential rotation until <j8] 

2.3. Energetics, and implications for the lightcurve and 

spectrum 

The stellar tidal acceleration creates regions of high 
and low gravitational potential in the longitudinal direc- 
tion. Insolation heats the fluid, and this heat can be 
tapped to perform work against the stellar tidal field by 
causing pressure gradients which move fluid from regions 
of low to high potential. The ultimate source of energy 
for tidal heating driven by the thermal tide is then the 
stellar radiation field. The rate of work done by the ther- 
mal tide is 



A(TT) 
work 



-AM— I (-VU) ■ dl ~ -— iV< 
m J m 



TT) 



(G) 



where U is the tidal potential, a = 2{n — f2), m = 2, and 
the integral is taken along a fluid trajectory dl = R p d<f> 4>. 
Here 7V( TT ) ~ AM § (-VU) ■ dl is the thermal tide 
torque on the heated layer. In torque equilibrium, the 
thermal and gravitational tide torques balance giving 
^Y(GT) _ _7v (TT) . The heating rate due to the gravi- 
tational tide is then 

fL N (GT) = _jLjy(TT) 

m m 
showing explicitly that the work done on the atmosphere 
by the thermal tide is converted into heat by the gravi- 
tational tide in steady state. 

Since the work done on the atmosphere is ultimately 
powered by the stellar radiation, the maximum power 
input is given by the rate of absorption of stellar flux, 

2 
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where is the stellar luminosity and P orD = 2tt/u is 
the orbital period. We will show in f|5]that the efficiency 
of converting stellar flux to heat by gravitational tidal 
dissipation is as high as ~ f — 10%, and increases toward 
the star. 

Some of the heat deposited by insolation is converted 
into work, leading to less thermal energy immediately 
re-radiated back out into space. In the (unphysical) 
limit of 100% efficient conversion of heat to work, and no 
subsequent conversion of kinetic energy back into heat, 
the atmosphere would be far cooler than the equilibrium 
temperature, and the thermal emission from the planet 
would be solely due to the flux coming out from the core. 
In torque equilibrium, the work done on the atmosphere 
by the thermal tide is converted into heat by the gravita- 
tional tide. We presume the gravitational tide dissipation 
occurs in the convective core, and entropy is efficiently 
mixed throughout. The temperature profile as a function 
of depth, latitude and longitude may then differ markedly 
from the case with no tidal heating. 



First consider the change in the lightcurve for thermal 
radiation. For small flux from the core, and short ther- 
mal time in the absorbing layer, there will be a large 
day-night temperature difference. If, for example, 1% 
of the insolation is deposited in the core by the grav- 
itational tide, then the day-night temperature ratio is 
(Pcore/^*) 1/4 ~ 0.01 1 / 4 ~ 0.3, where F corc is the flux 
emerging from the convective core. For small F COIC , the 
only possibility of high temperature on the night side is 
redistribution of heat by zonal winds. Tidal dissipation 
generating a large -Fcore provides an alternative means 
to reduce the day-night temperature contrast (see figure 
1111 for the luminosity exiting the core versus planetary 
radius). 

A change in the vertical distribution of heat sources 
also affects the vertical temperature profile at the depths 
where the spectrum is formed. Depositing heat in the 
core rather that at the photosphere leads to a radia- 
tive flux more constant with depth near the photosphere. 
This alteration of the temperature profile may in princi- 
ple be imprinted on the spectrum. 

3. TEMPERATURE PROFILES AND QUADRUPOLE 
MOMENTS 

In this section we describe the model for time- 
dependent thermal forcing of hot Jupiter atmospheres 
due to asynchronous rotation and/or an eccentric orbit, 
as well as the the resultant gravitational forces on the 
perturbed atmosphere. The full nonlinear equations are 
presented, and then approximated as a time-independent 
background and harmonic perturbations. Analytic solu- 
tions are developed for the background and perturba- 
tions, and compared to the full nonlinear solutions. 

3.1. Equations and geometry 

We consider a planet and star of mass and radius 
(Mp,R p ) and (M*,i?*), respectively. The planet orbits 
with separation D(t) and true anomaly which we 
will treat as a nearly Keplerian orbit with semi-major 
axis a, eccentricity e and mean motion n = [G(M* + 
M p )/a 3 ] 1 / 2 . We consider an atmosphere in uniform ro- 
tation with rate fl, and the spin and orbital angular mo- 
mentum aligned. Differential rotation will be discussed 
in $5] We will work in a non-rotating coordinate system 
whose origin is at the center of the planet, and whose 
axes are fixed with respect to distant observers. Spheri- 
cal polar coordinates (9, 0) are used, where the star orbits 
at the equator with a colatitude 9 — tt/2 and longitude 
cf> = $(t), and a fixed point on the planet rotates with an- 
gular frequency cf> = Q. The derivative comoving with the 
planet is then d/dt — d/dt + fld/dcf). The cosine of the 
angle x between the vertical and the vector to the star is 
cos x = sin 9 cos(0 — $). The day side is over the angular 
range — tt/2 < <fi — <f> < tt/2. Let z be the altitude above 
some appropriate reference level and y(z) = dz'p(z') 
the mass column above that altitude. The atmosphere 
is treated as being both thin and in hydrostatic balance, 
hence the pressure P = gy, where the surface gravity is 
g = GM p /R 2 p . 

We use an approximate treatment of radiative trans- 
port in which the time-dependent insolation is treated 
as a specified heat source. We approximate the trans- 
fer of thermal radiation by flux-limited diffusion. For 
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the thermal radiation, we employ the solar composition 
"cond ensed" phase Rosseland opacities of lAllard et al.1 
(2001), which include the effect of grains in the equation 
of state, but ignores their opacity, as is appropriate if the 
grains have rained out to higher depth. There is a promi- 
nent dip in k along isobars centered around T = 2000 K. 
This feature will be apparent in our num erical results. 

We use the equation of state from ISaumon et al.l 
(1995), with 70% hydrogen and 30% helium by mass. 
At low density, this is an ideal gas equation of state 
P = pki,T / /im p , where p is the mass density, T is the 
temperature, m p is the proton mass, k\> is Boltzmann's 
constant, and /i ~ 2.4 is the mean molecular weight for 
a mixture of molecular hydrogen and helium. 

The temperature and flux are found by solution of the 
heat equation 



dT _dF 

and the equation for flux-limited diffusion 



F = 



16g sb T 3 A dT 
k dy 



(9) 



(10) 



Here C p ~ 7fcb/2/Ltra p is the specific heat per gram, and 
e is the heating rate per gram. Since column, or equiv- 
alently pressure, is assumed constant for fluid elements, 
the thermodynamic relation ds/C p = dT/T — V a< \dP/P 
is simplified to Tds — C p dT. We ignore vertical fluid 
motion relative to constant pressure surfaces, as well as 
horizontal motions relative to the mean rate Q. The 
quantity A is the flux limiter, which allows a smooth 
transition between the optically thick regime, A = 1/3, 
and the optically thin regime A — > 0. For convenience, we 
use th e limiter prescription of iLevermore &: Pomraningl 
(1981) 



A = 



R 



6 + 3R + R? 



(11) 



where R — (A/ n)\dh\T / dy\ is the ratio of photon mean 
free path to temperature scale height. The boundary 
conditions are (1) the radiation free streams at small op- 
tical depths, F ~ 4<r s bT 4 , and (2) F — > F cole at large 
y, where -Fcorc is the thermal flux that emerges from the 
convective core. 

We refer to direct solutions of eql9landll0las "nonlinear 
solutions", as temperature changes are not assumed to 
be small. The nonlinear solutions will be compared to 
solutions of the linearized equations which are valid only 
for small temperature changes. 

To numerically so lve eg 151 and [Tffl we adapted the code 
of lPiro et alJ 0005). We use second-order finite differ- 
ence in y and backward difference in t for stability. Typ- 
ically 128 y grid points were used, but selected results 
were checked at higher resolution and found to be accu- 
rate. The time step is set so that the average temperature 
at each grid point changes by < 10~ 3 over the time step. 
A background model, which we discuss below, is used 
as an initial condition, and the initial time taken to be 
noon (cf> — <&). Transients due to the initial conditions 
are observed to relax over a few forcing periods. 

3.2. Time- dependent insolation and Fourier components 



The time-dep endent heat ing rate per unit mass is given 
by Beer's Law (|Liou 112002ft 



e(y, 9, 4>, t) = k±F+(D) exp 



K*y 
cos\ 



e(cosx) (12) 



where the constant specifies the column (y ~ k^ 1 ) 
at which the stellar radiation is absorbed, and F+(D) = 
o~shT+ (R+/D) 2 is the bolometric stellar flux at the sub- 
solar point on the planet. We present numerical results 
only for a solar-like star with M* = M Q , = 5780 K 
and R* — R@. However, our analytic analysis allows for 
rescaling to stars of different mass and luminosity. The 
step function 0(cosx) enforces heating only on the day 
side. In radiative equilibrium, the thermal flux exiting 
the top of the atmosphere is 

F(y = 0, 9, <f>, t) = F COIe + / dy e(y, 9, 0, t) 



= F core + cosxF*6(cosx). (13) 

We take the upper limit of integration to be y — oo since 
y = k^ 1 is a shallow layer in the atmosphere, and the 
integral converges exponentially. 

The gravitational tide can only couple to fluid per- 
turbations with quadrupolar angular dependence and a 
specific harmonic time dependence. Therefore, we now 
expand the insolation in terms of a Fourier series in lon- 
gitude and time. Since cosx depends on only the com- 
bination of longitudes -0 = — $, we write 



e(y,9,4>,t)=K i! F*{D) ^ g„ 



sin#/ 



where the integral 

g m {s)-- 



2tt 



tt/2 



dtp e 



-imip — s sec -0 



(15) 



-it/2 



contains a contribution only from the day side. We com- 
pute the integrals g m (s) numerically. The m = term 
is the average over longitude. For eccentric orbits, us- 
ing F+(D) = F- k (a)(a/ D) 2 , the time-dependence can be 
expanded as a sum of harmonic terms using the Hansen 
coefficients defined in eq |A13l of appendix [Al 



D 



= £ Xl m {e)t 



ikn t 



(16) 



k= — oo 



The heating rate in eq !12l can then be written in a Fourier 
series in both longitude and time as 



e(y, 9, cj), t) = K+F+(a) X k 

rnk 



,im<p—ik. 



(47) 



The (to, fc)'th term in the series has a forcing frequency 
ffmfc = kn — mfl in the co-rotating frame and kn in the 
inertial frame. The (to, k) = (0, 0) term is the time and 
longitude averaged heating rate, which we denote as 



s(y,9) = F^(a)K ir Xg (e)g (K^y/ sin 9). 



(18) 



Since Xg m (e) = S m o(l — e 2 ) -1 / 2 , there are no time- 
independent (k — 0) non-axisymmctric (to ^= 0) heating 
terms. Axisymmetric (to = 0) time-dependent heating 
does occur for eccentric orbits, due to the radial force. 
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Fig. 2. — Temperature versus pressure over one diurnal forcing 
period for "deep" forcing at the equator (9 = 7r/2). Note the lag be- 
tween maximum forcing (at noon) and maximum temperature (be- 
tween noon and 3pm), and also the pronounced "bump" near the 
base of the heated layer at Pbase = <?/ K * = 10 5 ' 5 dyne cm -2 . This 
plot was produced using the nonlinear simulations for a circular or- 
bit with P or b = 4 days, super-synchronous rotation P sp i n = 3 days, 
k* = 10 -2 ' 5 cm 2 g _1 , g = 10 3 cm s -2 , F CO rc = 10 4 erg cm -2 s _1 
and a solar type star. Noon is the time of maximum heating. 

For the time-independent background, L dx go(x) = 
1/tt leads to a background surface flux 

/>00 

F(y = 0,6)=F col . c + dye(y,6) 
Jo 

= F COIC + -F4a)X™(e) sin0 (19) 
and thermal luminosity 

L = R 2 p J dQF(y = 0,(9) 

= AvR 2 p (p corc + ^F t (a)j (20) 

in thermal radiation from the planet. We will use eafTHl 
to define background models, on top of which we solve 
for perturbations. We define the time-dependent pertur- 
bation to the heating rate as 

6e{y,6,4>,t) = ^F i< {a) T Xl m (e)g m (4^) e im *-fZt) 

' V sin / 

For e = 0, Xf. m (0) = 8 m k enforces k = m, and the forcing 
frequency becomes a mm — m[n — f2), where 27r/|n — £2| 
is the diurnal forcing period. 

The values of g m (s) in the heating function decrease 
with increasing m. While the diurnal (m = ±1) terms 
dominate the observable flux variation, they do not con- 
tribute to the torque since the tidal potential has no 
dipole component. The quadrupole moments responsi- 
ble for spin and orbital evolution are dominated by the 
semi-diurnal (m = 0, ±2) components. 

3.3. Discussion of temperature profiles 

Figures [2] and [3] show temperature versus pressure 
over one diurnal forcing period for a planet in circu- 
lar orbit around a solar-type star. The orbital period 
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Fig. 3. — Same as figure [2l but with "shallow" forcing = 
10 -0 ' 5 cm 2 g — 1 . The base of the heated layer is now at Pbase = 
g/K* = 10 dyne cm -2 . The temperature perturbations diffuse 
significantly below the base of the layer directly heated. Here the 
forcing period, 2ir/a, is large in comparison to the thermal time 
t t ji at y < k^ 1 giving temperature maximum in phase with the 
maximum heating at noon. A thermal diffusion wave penetrates 
to large depths, and lags in phase with respect to the temperature 
at y < ■ This lag in phase leads to a torque that attempts to 
increase the spin of the planet. 

is P or h = 4 days and the spin period P S pin = 3 days, 
giving a diurnal forcing period of 12 days. Each curve 
is labeled by phase in the heating cycle, where noon 
means heating is maximum, etc. Figure [2] illustrates a 
heating function which extends "deep", down to a base 
pressure Pbase = g/^ — 10 5 5 dyne cm -2 , while fig- 
ure [3] shows "shallow" heating down to Pbase = g/K* = 
10 3 ' 5 dyne cm~ 2 . Naively one may have expected the 
temperature to change significantly only above Pbase- 
This is approximately true for the deep heating case. The 
shallow heating case, however, shows large temperature 
perturbations extending a factor of ~ 30 deeper in pres- 
sure than Pbase- Heating below Pbase is due to diffusion, 
rather than direct absorption of stellar radiation. We 
will show that there are two limits for the forcing period. 
In the "high frequency limit" , appropriate for planets at 
large separation from the star, the temperature pertur- 
bations are large only for depths P < Pbase- In the "low 
frequency limit", appropriate for close-in planets, tem- 
perature perturbations have time to diffuse down below 
Pbase, so that the star torques a layer deeper than that 
directly heated by the stellar radiation. 

The "bump" in the temperature profiles in figure [2] at 
P = 10 5 — 10 6 ' 2 dyne cm -2 is evidence of an inwardly 
propagating thermal diffusion wave. T his bump is also 
seen in temperature profiles in figure 5 oflShowma n et aLI 
(|2008f) . although the authors d o not identify this feature 
as such. IShowman et alj (|2008[ ) have a far more detailed 
solution method, including nonlinear three-dimensional 
hydrodynamics, angle and frequency dependent radia- 
tive transfer, and detailed chemical abundances and 
monochromatic opacities. That our solutions are in qual- 
itative and rough quantitative agreement in the region of 
interest lends confidence to our results. 

Another important comparison to the more rigorous 
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Fig. 4. — Thermal flux at the surface versus time for the cases 
shown in figures [2] and \3\ These simulations are at the equator 
(9 = 7r/2), flux is normalized by the average surface flux -F C ore + 
-F*(a)/4, and time is in units of diurnal forcing period (12 days). 
The solid black line is the outgoin g flu x if no net heat is absorbed 
or emitted by the atmosphere (eq|13|l. The blue long-dashed line 
is the shallow heating case (figure [3}. The red short-dashed line is 
the deep heating case (figure [2). 

results of Sho wman et al.l ((2008) concerns the lag in the 
temperature profile, shown in their figure 3. Although 
their simulation is of a synchronous planet, they generate 
a super-rotating equatorial jet which mimics the super- 
synchronous rotation in our results. Those authors in- 
deed find the temperature reaches a maximum eastward 
of the subsolar point; maximum temperature lags max- 
imum heating. This lag is especially prominent at large 
depth (their bottom panel). We regard this as strong 
evidence that, while inclusion of fluid motion may al- 
ter our results quantitatively, the qualitative result that 
temperatur e lags, due to therma l inertia, is in agree- 
ment with IShowman et al.l (|2008l ) . A crucial omission 
in their simulations is that the tidal force from the star 
is not included. We argue, based on the temperature lag 
seen in their figure 3, that the torque on this thermally- 
generated quadrupole will act to generate asynchronous 
spin. 

Figure 2] shows the thermal flux exiting the surface, 
in units of the average thermal flux, for the two cases 
shown in figures [2] and [3l The shallow heating case shows 
larger flux variations and smaller lag time in comparison 
to the deep heating case. In both cases, the maximum 
in thermal emission lags the maximum in heating (at 
noon) due to thermal inertia. The Spitzer lightcurves for 
HD 189733b show maximum thermal emission just before 
secondary eclipse. Using the geometry from figure [IJ 
maximum emission occurs just before secondary eclipse 
for a planet with super- synchronous rotation, in order 
that the temperature reaches maximum to the east of 
the subsolar point, after it has passed through noon. A 
similar effect would occu r for a super- synchron ous zonal 
wind, as in the model of IShowman et all (|2008l ). 

One technical point about flux-limited diffusion is that 
it produces temperature profiles which increase inward at 
small optical depth. The "greenhouse" case shown in fig- 
ure[2]is captured in our solutions. Recall that greenhouse 



Fig. 5. — Effective column £ (eq |24|l . converted into units of 
fractional mass of the planet, for the deep heating case (red, figure 
[2)| and the shallow heating case (blue, figure [3). The solid lines 
are the full numerical result. Fourier transform of the solid lines 
gives the diurnal component (dotted lines) and the semi-diurnal 
component (dashed lines). 

heating is due to <C k; stellar irradiation deposits 
heat at large optical depth, leading to a temperature 
increase by a factor ~ (k/k*) 1 / 4 above the skin tem- 
perature (F*/4cr s b) 1 ^ 4 . There has been recent interest 
in shallow heating, ^> k, due to TiO/VO absorption 
at P ~ mbar, producing an inwardly d ecreasing tem- 
perature profile, the "stra tosphere" case (|Hubenv et al.1 
I2003t iFortnev et afll2008fh To model this effect would 
require angle and frequency-dependent solution of the 
transfer equations, beyond the scope of this paper. We 
note, however, that below Phase, the stratosphere case is 
expected to have temperature ~ (i^/ergt,) 1 / 4 , similar to 
our results with flux-limited diffusion. Since the torque 
is applied well below Phase, we expect our results to be 
qualitatively, and perhaps quantitatively correct. 

3.4. Quadrupole moments 

Appendix [A] contains a derivation of the secular 
changes in e, a and Q due to a quadrupole moment. 
In addition, we derive the change in total spin plus or- 
bital energy due to this quadrupole moment. Here, we 
outline the method of computation for the thermal tide 
quadrupole moment, as well as review Darwin's theory 
of secular evolution for gravitational tides. 

The "thermal tide" force results from the mass mul- 
tipole moments induced in the atmosphere by time- 
dependent insolation. The multipole moments are de- 
fined in terms of the density field that is the solution to 
eqs. Inland [11)1 We have 

Q7J\t) = J d 3 x p(x,ty Y; m (6,4>). (22) 

The density p{x, t) is computed from the temperature 
profile T(x,t) using the ideal gas law and hydrostatic 
balance P — gy. Since r e p(x,t) is a real quantity, 
and Y e * m = (— l) m Yi^ m , the moments must satisfy 

(2&7 } ) * = Q^li" 1 )" 1 - As the static and azimuthally 



symmetric background model has zero quadrupole mo- 
ment, it can be subtracted off the integrand. We then 
make a change of integration variable from radius r to 
column dy = —pdr, and approximate r ~ R p , yielding 
the expression 



q(TT) 



(t) = 



■.R 2 p +e I dnY; mK 



Inspection of eq[23] shows that the tide couples to an 
"effective column" 



>:(6>. 6. t) \ 



dy 



p(x,t) - p(y,0) 



(24) 



Figure [5] shows numerical examples of the time- 
dependent column £ for the same simulations as in fig- 
ures[5]and[3] The effective column (solid lines) in general 
leads the heating function, a consequence of the tempera- 
ture lag and the constant pressure assumption. We have 
taken a Fourier transform of the time-series for S, and 
show the diurnal and semi-diurnal components, which 
also lead. The diurnal component dominates the tem- 
perature and density response, but produces no torque 
since the tidal force has no £ = 2, m = ±1 components. 
The semi-diurnal response is smaller, by a factor ~ 2 in 
the present example, but dominates the torque. 

The expressions for the secular change in the planet's 
orbit require the Fourier transform in time of the 
quadrupole moments, which we perform numerically 
from our simulation output. We define the desired trans- 
forms 



Q 



(TT) 
£mk 



ikn I 



m— — oo 



where 



Q 



(TT) 
£rnk 



2tt 



2-k/t. 



dt 



ikntniTT) 



J t III 



(*), 



(25) 



(26) 



and satisfy Q^Jl, _ h 



,(tt; 



For the non- 



(-i) m (q{ 

linear simulations, we evolved the atmosphere for 10 forc- 
ing periods, and computed the Fourier transforms numer- 
ically using orbits 9 and 10. 

We now turn to the gravitational tide. The stellar grav- 
itational tidal acceleration — WU induces fluid flow in the 
planet. Dissipation causes a phase lag ~ l/Q between 
the acceleration and high tide, where Q is the tidal qual- 
ity factor. This phase lag is the origin of the "gravita- 
tional tide torque" , wh i ch att empts to make the planet's 
spin synchronous. iHutl (|1981[ ) contains a clear review of 
Darwin's theory of tidal evolution. There it is assumed 
that the fluid response lags forcing by a time ri ag . The 
quadrupolar fluid response can be viewed as two point 
masses each of mass (k p /2)M i ,[R p /D(t— ri ag )] 3 , with lon- 
gitude 4>(f - ri a g) - - n ag ) and <&(f - r lag ) - f2(i - 
T ia g ) + TT- Here k p is the apsidal motion constant of the 
planet, which takes into account the internal structure 
in the external potential perturbation produced by the 
quadrupole moment. Under the assumption that ri ag is 
small in comparison to the forcing period, we expand to 
first order in ri ag . We ignore the term independent of 



Ti a g as it is in phase with the stellar gravitational tidal 
forcing and therefore, cannot lead to a torque. We find 
the effective gravitational tidal quadrupole moment 



Q'£?\t)=iT lag k p ALR 2 p Y 2m (Tr/2 : 0) ( 5e 



Y^{kn-mtt)X 2 k m {e)e 



ikn I 



(27) 



which is defined in the inertial frame. Rewriting 
this expression in terms of Q' p = (3/2fc p )(l/nri ag ) 
(|Mardling fe Linll2002D and by taking the Fourier trans- 
form in time, we find the imaginary component 




xY 2m (7r/2,0)X z k 

As a chec k of this result, we compared eq|28l IA161 IA17I 
and IA181 against Hut's analytic formulas for a range of 
orbital period, spin period and eccentricity. 

Note that Hut's "constant time lag" prescription for 
orbit and spin evolu tion disagrees with t he "c onstant lag 
angle" approach of iGoldreich fc Soterl (JT966). The use 
of a constant lag angle, with a constant Q, is unphysical 
as the torque is discontinuous as the spin changes from 
sub- to super-synchronous and vice- ver sa. By adopting 
the lag in time approach of lHutl (|1981h or lMardling fc Lira 
(2002) with the ab ove choice for Q, we find agreement 
with the results of IGoldreich fc Soterl (|1966l ) when the 
planet is far from a synchronous spin state, as in the 
case of the Jupite r-Io system. There f ore, th e constraints 
on Q inferred by IGoldreich fc Soterl (Il966f) may still be 
applie d. N ote, however, that llngersoll fc Dobrovolskisl 



1980T 1 utilize the 



(|1978h and iDobrovolskis fc Ingersoll 
constant lag angle approach. We find both qualitative 
and quantitative differences with their work. 

3.5. Linear Perturbation Theory 

For large forcing frequencies, the time-dependent tem- 
perature and flux changes are small compared to the 
time-average values. One may then treat the time- 
dependent changes 8T(y,0,<fr,t) and 8F(y,0,<f>,t) as lin- 
ear perturbations about the time-independent back- 
ground values T(y, 9) and F(y, 9). 

The background model is computed by integrating 



dT 
dy 



nF 



16a sb T 3 A 
with flux found by integrating eg J18I 

F(y 7 0)=F COIC + dy'e(y',6) 



(29) 



= F a 



>cos0 e 



and subject to the boundary condition T ~ 
(F(Q, 9)/4a s b) 1 / 4 at small optical depth. We construct 
this model by integrating inward from the surface. An 
example is given in figure [S] Our treatment of radia- 
tive transfer is simplified, but nevertheless reproduces the 
main features of full solutions to the transfer equations. 
At small optical depth the solution becomes isothermal 
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Background model constructed using time- 



usmg re* 

10 4 erg cm -2 s , g = 10 3 cm s~ 2 , a solar-type star and P OJ 



10" 



1 

* 3 



at the "skin" temperature. Heating by inward-going 
stellar photons, here through the heating function, gen- 
erates an equal outward flux of thermal photons (for 
the time- independent case). At a column y ~ k 
the heating function and outward flux decrease exponen 
tially, and the temperature profile becomes isothermal at 
T^- 1 , 9) ~ (F(0,60/o- sb ) 1/4 (K/K*) 1/4 - Lastly, the flux 
-Fcorc from deep in the atmosphere eventually causes the 
temperature to increase inward, the gradient becoming 
steep enough to cause convection. 

The linearized equations are found by perturbing eq|5] 
and [10] while enforcing strict hydrostatic balance for 
the perturbations, so that each fluid element remains 
at constant y. We eliminate dependence on (f> and t 
by expanding ST and SF in in a Fourier series (e.g., 
cqll7j) . The time dependence of the perturbations is then 
d/dt = —ia rn k- We will suppress m and k in the following 
equations. The linearized form of eg 151 and 1101 are 



dST _ dT/dy 
dy 



1-A, 



SF r A / 

— + {k t - 3 - A r (k t 



^ ST 

1)}jT< 



(31) 



and 



dSF 
dy 



= —iaCpST — Se, 



where kt = <91nK/<91nT|p and 

i? 2 (i? + 4) 



A, 



(6 + 3i? + R 2 ){2 + R)' 



(32) 



(33) 



At large (small) optical depth, — * (1). We solve 
eqEUand EH subject to the boundary conditions 5F/F = 
A5T/T at the top of the grid and 8T,6F -> at the 
base. The equations are solved by finite difference in 
y to obtain a banded matrix equation which is readily 
inverted to find the response ST and SF to the heating Se. 
The top of the grid is set at a column 10 -2 times smaller 
than the thermal photosphere, and the base of the grid 
is set at the radiative-convective boundary, typically well 
below y — n^ 1 . 



Fig. 7. — Semi-diurnal (m 
bations for a planet at P or b 
rate of P„ 



= 2), linear temperature pertur- 
4 days with a prograde rotation 
spm — 1-0 day, giving a semi-diurnal forcing period of 
2/3 day. Surface gravity is g = 10 3 cm s — 2 , 8 = 7r/2, and 
= 10 -2 5 cm 2 g _1 . The high frequency limit from cq |38l is 
a good approximation for P > 10 4 dyne cm -2 . 
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Fig. 8. — Same as figure [7] but with a spin period P sp j n = 
3 days and shallow heat deposition re* = 10 — 5 cm 2 g _1 . The low 
frequency limit is a good approximation for P < 10 6 dyne cm -2 , 
below which the oscillations in the real and imaginary part are 
evidence of a diffusion wave. 

Given the temperature perturbation, the density per- 
turbation is Spjp = —ST/T at constant pressure. To 
linear order, the effective column is given by 



E(0): 



dy 



ST{y,( 
' T(y,6 



(34) 



implying that the temperature perturbations lagging the 
heating tend to torque the planet away from synchronous 
rotation. 

Examples of the solution of the linearized equations are 
shown in figures [7| and [HJ along with analytic approxima- 
tions to be discussed presently. The difference between 
these two plots is mainly due to the depth of heating, i.e. 
deep in figure [7] and shallow in figure [8] 
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For the deep heating case in figure[7l the forcing period 
(=2/3 day) was chosen to be shorter than the diffusion 
time for the base of the heated layer. There is no time 
for diffusion to transport heat, the dF/dy term in eqj9]is 
small, and the temperature fluctuates solely due to the 
time-dependent heating. A longer forcing period leads 
to larger input of heat, and effective column E a Pj. 
where Pf = 2it/\(t\ is the forcing period. An analytic 
solution will be presented in £)3.6l Since the imaginary 
part of the temperature perturbation is large only above 
-Pbase = the torques will extend down to Pbasc- 

The shallow case in figure [8j on the other hand, 
has a sufficiently long forcing period (6 days) that the 
perturbations extend deeper than Pbase by some 2 or- 
ders of magnitude in pressure. For pressures P < 
10 6 dyne cm -2 , the first term in eq[S]can be ignored due 
to the small forcing frequency. In this limit, the atmo- 
sphere responds in phase to the forcing, changing quasi- 
statically from one equilibrium to another. No torque is 
exerted on the atmosphere in this limit since the temper- 
ature perturbation is in phase with the heating function. 
That is, for the semi-diurnal mode, temperature maxima 
occur at noon and midnight, while temperature minima 
take place at sunrise and sunset. In this "cr = 0" limit, 
we can solve for the (real) temperature perturbation by 
integrating eql9]with dT/dt — 0. In this regime, the ther- 
mal inertia of the material is small is comparison to the 
amount of photon energy absorbed per cycle and there- 
fore, the temperature is "locked" to the heating func- 
tion. For the forcing frequencies of interest, this limit 
will always be invalid sufficiently deep in the envelope, 
at which point the solution is well described by a diffu- 
sion wave. This wave is apparent in figure [8] in the region 
P = 10 6 — 10 7 dyne cm~ 2 , where the real and imaginary 
parts are comparable, and exhibit oscillations. The first 
and second terms in eq. [9] are comparable for diffusion 
waves below the base of the heating layer. The decay of 
the envelope of this wave into the planet determines the 
column to which torque is applied. 

Using dimensional analysis on eq|9] and I10[ heat can 
diffuse down to a column 



Udm : 



3kCJ<t| 



(35) 



in a forcing period Pf — 2%/\a\. This depth is to be com- 
pared to the base of the heating layer at k^ 1 . Heating 
can diffuse below the base of the heating layer when 

> _ 3tt kC„ 



= 3.6 days 



K 



0.01 cnr g 



0.1 cm 2 g 
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2000 K 



k* I \ T 
In the low frequency limit, the effective column 



(36) 



,1/2 



(37) 

than the high 



£~2/diff oc P f ' 

which increases more slowly as Pf 
frequency limit. 

Close to the star, where the gravitational tide is large, 
the equilibrium spin must be in the long forcing pe- 
riod limit, and far from the star, where the gravitational 



torque is small, the equilibrium spin must be in the short 
forcing period limit. The dividing line between these two 
limits depends on factors such as stellar flux, opacities 
and tidal Q. 

Unless dT/dy = T/y, the "diffusion depth", ydiff, 
differs from the often-used "cooling depth" , y C ooi = 
F/(C p T\a\). See the lower panel in figure El for a nu- 
merical example. In the problem at hand, y coo i and j/diff 
can differ by orders of magnitude, especially at small 
optical depths and below the base of the heating layer 
where the flux drops exponentially. Physically, y coo i is 
the depth down to which the heat content can be ra- 
diated in the timescale 2tt/|ct|. In the absence of heat- 
ing, the flux t ends to become constant with depth for 
V % J/cooi (e.g. iPiro et al"1l2005f h As we will show in the 
next sections, j/diff is the depth down to which small per- 
turbations can propagate, and hence is the more relevant 
for heating deep in the atmosphere. 

3.6. Analytic solution in the high frequency limit: 
o *th > 1 

The thermal time at the base of the absorbing layer, 
£th ~ C p T/ K+F+, is the timescale over which the layer 
can heat or cool. For large forcing frequencies atth ^S> 1, 
heat cannot diffuse over significant distances in a forcing 
period and the term ddF/dy can be ignored. In this limit, 
the temperature perturbation is determined solely by the 
local heating function (eq l32"|) . i.e., 

STmk^-^r, (38) 

where 5e m k is the appropriate Fourier coefficient in eq |171 
This expression shows that the temperature perturba- 
tions are small in the high frequency limit, ST decreases 
exponentially below y — k" 1 and its magnitude in- 
creases with forcing period since there is more time to 
absorb heat. With respect to the forcing frequency a m k, 
the temperature perturbation lags the forcing by 90° in 
phase. It is compared to the solution of the linearized 
boundary value problem in figure [7] 
The effective column is then 

' 5e mk {y,0y 



£ mit (0)~ — 



0V, 



dy 



T(y, 



(39) 



Since T(y, 9) varies relatively slowly for y < ft" 1 , we ap- 
proximate it as a constant and pull it out of the integral. 
The remaining integral is 



ds g m (s) = G r , 



(40) 



where Go = 1/tt and G2 = 1/3tt. Plugging eqflOl into 
eqj39]we find 



o~ 



C n T 



mk^p 



For the diurnal (semi-diurnal) component, eg 14 II implies 
£ leads the heating function by 7r/2 (tt/4) in longitude. 
This phase lead corresponds to 6am (9am) for the max- 
ima of the diurnal (semi-diurnal) components in figure 



The quadrupolc moments require the angular integrals 



F lm = 2n / d6sm 2 eY lm {0 : {)). 
Jo 



(42) 
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For example, F22 = \/l57r/2(37r/16) and F 2 q = 

(— 7r/16)v57r- The final result for the quadrupole mo- 
ments in the high frequency limit is 



stellar fluxj 



R 2 ^ 



tmk J 



iCpT 



Xl m (e)F £m G m , (43) 



confirming the order of magnitude estimates found in 
eq[3] and |1 if we set AM = R 2 / 'k*, i t h = C p T / K+F^a) 
and 5=1. 

Eq|43 | agrees with the quadr upole moment found in 
llngersoll fc Dobrovolskisl |l978), up to the uncertainty 
as to what value of temperature should be used when it 
is pulled out of the integral. These authors derived the 
high frequency limit appropriate for a circular orbit and 
asynchronous rotation, and implicitly assumed a finite 
atmosphere bounded below by a hard surface. Eql43l 
generalizes their result to a deep atmosphere, and allows 
for an eccentric orbit. 

3.7. Analytic solution in the low frequency limit: 
(J t t h < 1 

We now consider the low frequency limit atth <C 1, 
in which the thermal time in the layer directly heated 
by insolation (y < k^ 1 ) is short compared to the forc- 
ing period. The small thermal inertia of the absorbing 
layer implies that temperature perturbations are large 
and in phase with the forcing, with corresponding den- 
sity perturbations 180° out of phase with the insolation, 
implying that the torque is zero. There will always be a 
layer at y > ft* 1 sufficiently deep that the thermal time 
becomes longer than the forcing period, and the temper- 
ature perturbation in this deep layer will be out of phase, 
yielding a torque. 

To illustrate this, we divide the atmosphere into two 
separate plane-parallel layers, as depicted in Figure M 
The stellar flux is absorbed in the upper layer at a col- 
umn y ~ k^ 1 , and generates an outward flux of thermal 
radiation. In the lower layer, the stellar flux is negligible, 
and time-dependent temperature perturbations can only 
penetrate this layer due to thermal diffusion. 

The thermal time in the upper layer is short, and the 
time derivative term in eq. [3H can be set to zero. The 
temperature profile changes quasi-statically in response 
to changes in the heating rate. As an analytic example, 
consider a flux profile F(y) — Fq exp(— K+y) and con- 
stant opacities k and ft* -C ft. The temperature deep in 
the layer is then T ~ (3fti 7 b/4<7 s bft*) 1/ ' 4 . If the heating 
produces a change in flux Fo — > Fn + SFq , the fractional 
change in temperature deep in the upper layer is 



T 1 )- 4iV 



(44) 



For semi-diurnal forcing, the fact that ST reaches it's 
maximum at noon and midnight implies that 5p, or the 
thermal tidal bulge, achieves it's maximum value at sun- 
rise and sunset. As result, the gravitational tidal field 
cannot apply a torque to the the thermal tidal bulge in- 
duced in the upper layer in the o~t t h <C 1 limit. 

Direct heating is negligible in the lower layer, and Se 
in eq. [321 may be set to zero. Therefore, the tempera- 
ture perturbation is determined by a source- free diffusion 
equation subject to the boundary condition in ea !44l at 
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Fig. 9. — Cartoon depicting time-dependent heating in the 
limit where the forcing period, 2ir/a, is large in comparison to the 
thermal time t t ^ of the absorption layer. The temperature of the 
absorption layer is approximately "locked" to the time-varying in- 
solation such that the maxima in temperature take place at noon 
for diurnal heating and noon and midnight for semi-diurnal heat- 
ing. However, a thermal diffusion wave penetrates to a column 
depth y ~ j/diff = k-* 1 /V a *th- As expected for a thermal diffusion 
wave, the temperature fluctuation at y ~ J/diff ia g s t ne temper- 
ature fluctuation in the absorption layer by n/A in phase, which 
ultimately allows the heated layers to couple to the semi-diurnal 
tidal field in such a way as to increase its spin angular momentum. 

y = « ~ 1 . The background temperature T is roughly con- 
stant below y = ft* \ and we assume constant opacity ft 
for simplicity. With this, the temperature perturbation 
in the lower diffusion layer obeys 



i*C<6T^ 6 ^ T3d2ST 



3 ft 



dy 2 



(45) 



over the interval ft* 1 < y < 00. The general solution of 
eql45l has the form 



where 



-i7r/4 



16(7sbT 3 



(47) 



is the wavenumber, in units of inverse column. Note that 
a may be either positive or negative in this expression. 
Finitcness of the solution at y — 00 implies C = for 
either sign of a. The final result for the temperature 
perturbation in the lower layer is then 



ST. SF _ 
-T {V) = AF E 



s(y- 



) 



(48) 



Since s is complex, eq 1481 implies oscillations in tempera- 
ture, as well as an envelope decreasing to larger depths. 
These oscillations are apparent in the linear solution in 
figure [H as well as in the nonlinear solution in figure [H 

The effective column is determined by integrating eq|48l 
over the lower layer to find 



£ = - 




(49) 



which leads in phase by 37r/4, as compared to tt/2 for 
the high frequency limit. For the diurnal (semi-diurnal) 
component, eq |49l implies S leads the heating function by 
37r/4 (37r/8) in longitude. This phase lead corresponds 



12 



to 3am (7:30am) for the maxima of the diurnal (semi- 
diurnal) components in figure O 

Combining the results from H3.6I and this section, we 
find generically that the semi-diurnal component of E, 
and hence the quadrupolc moments, have the opposite 
phase as the gravitational tide quadrupole moments. 
Hence they promote asynchronous spin, and may drive 
eccentricity. Furthermore, cq]49l shows that £ continues 
to increase as the forcing frequency goes to zero, albeit 
with a shallower power of a. We summarize the scalings 
from cq[4T] and [49] by writing 

E ~ k- 1 Min [(a t t h) -1 , {a ith)" 1/2 ] , (50) 

and phase lead in the range 7r/2 — 37r/4. 

4. EQUILIBRIUM SPIN RATE 

Hot Jupiters are commonly assumed to be very nearly 
synchronized. The argument is that the synchronization 
timescale, for Q' p = 10 6 , is far shorter than the age of 
observed planets. The torque on the planet due to the 
gravitational tide is given by (eq[28] and IA18[) 




ft) 



n 



~4i 



-7 x 10 



for a circular orbit. Starting from a rapid initial spin rate 
Oq n, the synchronization time is 
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(52) 



where I p = r/M p R p , rj ~ 0.25 is appropriate for Jupiter- 
like planets, and we have defined the "dynamical time" 
P p = 2tt{RI/GM p ) 1 / 2 . Hence, for Q' p = 10 6 , planets out 
to P or b ~ 150 days should be synchronized in 5 Gyr. 

The fact that gas giants have Q' p ~ 10 6 implies gravita- 
tional tides are far weaker than for terrestrial type plan- 
ets with Q' p ~ 10—100. Additional torques are relatively 
more important for gas giants in comparison to terres- 
trial planets. We now show that thermal tidal torques 
are capable of generating significant asynchronous spin 
for hot Jupiters. 

The high frequency limit has a simple analytic solu- 
tion. For zero eccentricity, plugging the £ = m = k = 2 
quadrupole moments from eq ]43l into eq |A23l gives the 
thermal tide torque in the high frequency limit 



at( tt ) = 



3tt 



n 2 B*F* 



16 (fi -n)C p T 



(53) 



Recall that T should be evaluated at y ~ k+ 1 , and 
here F+ = F+(a). This estim ate agrees with that of 
llngersoll fc Dobrovolskisl (|1978h . Plug ging in numbers 
for a solar-type star we find 



TT) 



= 2.8 x 10 32 



erg 



4 days 



where R w = R/10 10 cm and f T = T '/ '(F*/ 'cr sb ) 1/4 . 
Clearly for the hot Jupiters the gravitational tide torque 



in eg 1511 and the thermal tide torque in eg 1541 are com- 
parable. Equating eq !54l and 1511 we find the equilibrium 
spin rate 



O e q = U 1 T 



In F* GQ' p 
24 C P T n 3 R p 



(55) 



or, numerically 



Pn, 



1± 



0-63 (fe) (^R^fr 1 



T7n (56) 



The =F sign denotes the subsynchronous (Q < n) and su- 
persynchronous (fl > n) solutions. The subsynchronous 
solution becomes retrograde (CI < 0) for sufficiently long 
P or b- Eql56l shows that the degree of asynchronous spin 
increases for larger P or b and is large for the hot Jupiters 
for Q ~ 10 5 — 10 6 . At long orbital period, the spin period 
asymptotes to a constant 



Pm 



±6.34 days ( R 



p f-. 

10 1q6 Jt 



1/2 



(57) 



If no new physics intervenes, this asymptotic spin rate 
should hold for planets at P or b ^150 days which have 
sufficient time to attain spin equilibrium. Since the de- 
gree of asynchronous spin decreases toward the star, the 
high frequency limit becomes less accurate there. 

Small perturbations about the equilibrium spin rate, 
at fixed P or b and R p , are stable for both the super- and 
sub-synchronous solutions given above. 

The accuracy of the high frequency limit can be esti- 
mated analytically. The semi-diurnal (m = 2) forcing 
period in the high frequency limit is 



P/ = 3.2 days R 



1/2 
10 



10 e 



6\ V2 



(58) 



independent of orbital period. The diffusion time to the 
base of the heating layer is (eqC 



P diff ~1.6 days 



orb 



4 days / 



10" 2 - 5 cm 2 



.(59) 



At long orbital periods, the forcing period becomes 
shorter than the diffusion time, and the high frequency 
limit is applicable. Equating the expressions in ea. 1581 
and [59] a rough estimate for the applicability of the high 
frequency limit is 



Porb > 8 days P^ 2 



10 e 



6\ 1/2 



(60) 



The high frequency limit allows useful analytic intu- 
ition, but neglects diffusion and approximates temper- 
ature at the base as constant. We now present results 
on equilibrium spin rate and heating rate relaxing these 
assumptions. Note, however, that these calculations are 
not fully self-consistent in that we fix a radius, and we 
also take, for simplicity, the flux from the core to be 
constant, independent of latitude. In SJ7] we will self- 
consistently determine the radius in thermal equilibrium 
and determine the self-consistent value of P CO rc- 

Figure [Tol shows equilibrium spin rate for three differ- 
ent methods of calculation: "nonlinear simulations" as 
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Fig. 10. — Equilibrium spin period as a function of orbital 
period for super-synchronous solutions (f2 > n). The "shallow" 
and "deep" lines represent ft* = 10 -0 ' 5 , 10 — 2 5 cm 2 g _1 , respec- 
tively. The solid black line is the analytic model from ij3.6l with 
fx = 1. The short-dashed green line is the orbital period, and 
the long-dashed green line shows P sp m = (ll/18)P or b> the spin 
period below which the gravitational tide excites eccentricity. Sur- 
face gravity is g = 10 3 cm s~ 2 and radius is R v = 10 10 cm. The 
flux from the core is taken to be -F C ore = 10 4 erg cm -2 s — 1 , inde- 
pendent of latitude. 

discussed in £13. II and 13.31 "linear" calculations solving 
the linearized equations as a boundary value problem 
( §3.5|) . and the "analytic" result from eql56]in the high 
frequency limit. Two different heating depths are shown, 
"shallow" and "deep" lines, as well as two different values 
for the tidal Q. The good agreement between linear and 
nonlinear calculations implies that at the large depths 
where the quadrupole moments are determined, nonlin- 
ear effects do not greatly change the result. The analytic 
approximation is accurate to with ~ 10%. The degree 
of asynchronous spin is large for both Q' p = 10 5 and 10 6 
cases. The shallow heating case is more nearly synchro- 
nized than the deep heating case as it is more in the low 
frequency limit than the deep heating case, and hence 
the thermal tide torque is weaker. 

5. TIDAL HEATING RATES 

The large tidal heating rates for eccentric hot Jupiters 
can inflate the radii of these planets if heat is deposited 
in the convective core, or slow thei r cooling if heat is 
depo sited in the radiative envelope (Bodcnhe imer et al.l 
[20011) . Previous investigations considered transient tidal 
heating due to an energy reservoir in the initial eccentric- 
ity, or steady-state heating due to eccentricity pumped 
by a companion. In this section we compute gravitational 
tidal heating rates due to the asynchronous spin set by 
a balance of gravitational and thermal tide torques for 
a planet on a circular orbit. We compare these heating 
rates to the cooling rate of the core to determine if tidal 
heating can power the large observed radii. 

Appendix [B] describes an updated ve rsion of the plane- 
tary s tructure and evolution models of lArras fc Bildstenl 
( 2006) . Figure QT] shows the cooling luminosity of the 
convective core as a function of radius for P rb = 2 days 
around a Sun-like star. Also shown is the timescale to 
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Fig. 11. — Cooling luminosity from the core (curves sloping up 
to the left) and cooling time (curves sloping up to the right) versus 
planetary radius for three different planetary masses, M p /Mj = 
0.5, 1.0, 2.0. A solar type star is used and the planets sit at P or b = 2 
days. The heating depth is given by ft* = 10 — 2 cm 2 g _1 . Closer 
to the star, insolation is larger and L COTC decreases, and visa versa 
further from the star. 

attain thermal equilibrium 



Pcorc 

Cooling timescales for Jupiter-mass planets become 
shorter than a typical age of 3 Gyr for R p > (1.3 — 
1.4) x Rj. Cooling rates closer (further) from the star 
are lower (higher) due to the insulation effect. 

Models for the passive cooling of irradiated hot 
Jupiters can typically explain planets with radii R < 
1.2Rj and below, but have difficulty slowing the cool- 
ing enough t o explain planets with R > 1.2Rj (e.g. 
lFortnevi r2008). If we are interested in the radius range 
R = (1.2 — 1.8) x Rj, figure [TT1 shows this corresponds 
to core cooling rates of 10 25 — 10 28 erg s _1 , depending 
on mass and insolation level. In order for tides to halt 
the cooling of the core, this loss of heat at the radiative- 
convective boundary must be balanced by a gain of heat 
in the core due to tides. 

The tidal heating rate in eq |A21l is related to the torque 
in eq |A18l by the pattern speed (kn — mQ,) Jm for each 
harmonic. In spin equilibrium, with _/V ( ' TT ) = — N^ GT \ 
the gravitational tide heating rate for the high frequency 
limit is 

= 2.9 xlO 27 ergs" 1 < (62) 

\ Porb / 

The tidal heating rate in the high frequency limit is inde- 
pendent of Q since fl—n cx (Q'p) 1 / 2 and E cx (Q— n) 2 /Q' p . 

Figure [12] shows the tidal heating rates for the equilib- 
rium spin rates in figure [TU1 The same small deviations 
between the linear and nonlinear results are apparent in 
this figure, but note that for the deep heating case the 
nonlinear, linear and analytic (high frequency limit) re- 
sults agree at large orbital period. The shallow heating 
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Fig. 12. — Gravitational tid e heating rate for the equilibrium 
spin periods shown in figure ITUl 

case shows a smaller heating rate than the deep heat- 
ing case as the torque is smaller. Recall that the tidal 
heating rate is independent of Q' p in the high frequency 
limit. This agrees with the deep heating case for large 
Porb- Somewhat counterintuitively, larger Q' p leads to a 
larger tidal heating rate at smaller P or b when thermal dif- 
fusion effects become important. The dependence of spin 
and heating rates on Q' p can be understood simply. Let 
AT(TT) = a ( n )(fi- n )- a and iV( GT ) = -b{n){tt-nfQ-\ 
for some functions a and b of orbital period, and expo- 
nents a and (3. The equilibrium spin rate is 



a. 



aQ 
~b 



l/(a+/3) 



(63) 



and the tidal heating rate is 



E G = (n- n y +p bQ- L = a 



^PhD- 1 = „(/3+ 1 )/(/ 3 + Q ) 



(Q/b) 



(l-a)/(/3 



In the high frequency limit, for Darwin's theory of tides, 



1 and (i ~ 1 give f2 e 



Q 1 / 2 and E oc Q°. 



In the low frequency limit, diffusion leads to a ~ 1/2 
so that fi e q — n oc Q 2 ^ 3 and E oc Q 1 ^, i.e. heating 
rate increases with Q. Lastly note that even varying 
Q by an order of magnitude and the heating depth by 
two orders of magnitude the dissipation rate varies by 
less than one order of magnitude, and is in the range 
E ~ 10 26 -10 28 erg s" 1 for hot Jupiters at P orb < 5 days. 
This heating rate is sufficient to explain the large radii 
of hot Jupiters, as seen by comparing figures fTTI and [T"2l 
Thermal equilibria are stable to small perturbations 
in R p if the perturbation to cooling dominates that of 

heating. If P hoat oc Rp and L corc oc R p , thermal stability 
requires d > c. The tidal heating rate in eqjBUhas c = 4. 
Figure [TT] shows that the cooling rate has d ~ 8—16 over 
the range of radii shown, implying stability. We find, 
however, that at a larger radius R p > 3Pj, the exponent 
d becomes less than 4, implying thermal instability. 

6. ECCENTRICITY EVOLUTION 

The orbits of some hot Jupiters are clearly more circu- 
lar than for the population of long period planets. This 



is often attributed to circularization by the gravitational 
tide raised in the planet by the star. Th e rate of change 
of eccentricity is (jMardling & Linll2002t) 
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For the equilibrium spin rate implied by N^ GT ^ — (see 
Hut (1981), eq.ll), the bracketed factor in eq[65]is pos- 
itive, implying circularization. The circularization time 
is then 
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4 hrs 
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V103MJ {Me)~^Me) 
For hot Jupiters with radii R = (1.0— 1.8) xPj 



10 6 



(68) 



< 10 6 



and P or b < 5 days the orbits should be highly circular 
after 5 Gyr. While many orbits are co nsistent with zero 
eccent ricity, it has been emphasized bv lMatsumura et al.l 
( 2008) that nearly 1/4 of planets within 0.1AU have e > 
0.1. In the absence of additional perturbations to the 
orbit, the finite eccentricity can only be explained by 
rather large Q' p ~ 10 s . The question is then why planets 
with presumably similar structure and orbital periods 
should have tidal Q differing by orders of magnitude. 

Inclusion of both the thermal and gravitational tides 
changes this simple picture in two ways. 

First, eqEn] shows that there is a critical spin rate 
O = n(18/ll)(/i//2) above which the gravitational tide 
excites, rather than damps, eccentricity. For evolution 
only under the gravitational tide, such rapid spin rates 
occur only for young, rapidly rotating planets. Synchro- 
nization occurs on a timescale shorter than the circular- 
ization time by ~ 10 3 , so that negligible eccentricity is 
excited in the early spin down phase. However, taking 
into account both gravitational and thermal tide torques, 
2]shows that significant deviations from synchronous ro- 
tation occur in the steady state, over the entire life of the 
system. For e = spin equilibrium, figure fTOl shows the 
critical spin frequency ft = (18/ll)n. For Q p = 10 6 , 
we find the gravitational tide can drive eccentricity for 
Porb ^ (5 — 9) days, depending on the depth of the heat- 
ing layer. For Q' p = 10 5 the critical orbital period is 
longer by a factor of a few. 

Second, the thermal tide may excite eccentricity due 
to the fact that density leads the forcing. The small 
eccentricity (e -C 1), high frequency (atth 3> 1) limit is 
analytically tractable. To le ading order in e, the results 
of ^3.61 plugged into eq |A24l yield 
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Fig. 13. — Growth rate of eccentricity e in the e <C 1 l imit, 
using the high frequency limit for the thermal tide (scccq |70t and 
the equilibrium spin frequencies for e = from eq |56l Positive 
values denote eccentricity growth and negative values denote ec- 
centricity damping. The solid (dashed) lines are the super- (sub-) 
synchronous solutions for the spin. The upper panel is for strong 
gravitational tide (Q' p = 10 5 ) while the lower panel is for weak 
gravitational tide (Q' p = 10 6 ). A solar type star was used, and 
M p = Mj and R p = X.3Rj. The divergences are due to the 
(fc , m) = (3, 2) and (1, 2) resonances for the thermal tide terms in 
eq|70l 

The four terms in the last parenthesis are the 0(e) 
harmonics (m,k) = (0,1), (2, 3), (2,1) and (2,2). De- 
pending on fi/n, eccentricity can either be excited or 
damped. Large response occur near the three resonances 

Pspin/ -^orb = 2/3, 1, 2. 

At small eccentricity, both the gravitational and ther- 
mal tides have e oc e, implying either exponential growth 
or decay with the rate (see eqESl HH1 and 155]) 
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Figure [TBI evaluates eq[70] using £1 = fl eq for e = from 
cq 1561 for the supersynchronous solution (solid line) and 
the subsynchronous solution (dashed line) . Two different 
tidal Q' p = 10 5 - 10 6 are shown. The value of Q' p deter- 
mines how long the orbital period must be before the spin 
becomes sufficiently asynchronous for 7 e to change sign. 
In both cases shown in figure [T3l the same basic pattern 
is found. For the supersynchronous case, close to the star 
eccentricity is strongly damped, then a region of eccen- 
tricity growth ending with a resonance, then a region of 
eccentricity damping starting with the same resonance. 
The eccentricity is always driven at long orbital periods, 
although the growth rate is too small to affect the orbit 
over the system's lifetime. In the regions of large ec- 
centricity growth, even if a planet began with small e, 
10-100 cfoldings would imply that the eccentricity could 
grow to large values. For small Q' p = 10 5 , eccentricity 
may be driven even out to P Q rb = 10 days, implying the 



observed eccentricity at these relatively large orbital pe- 
riods may not be primordial, but rather influenced by the 
thermal tide. For the subsynchronous case, the regions 
of eccentricity growth are far more limited, and may be 
confined to orbital periods very near the resonance. 

We have shown that small, but finite, eccentricity is 
unstable for small "windows" in P OI b- Hence even if all 
orbits began as nearly circular, a non-monotonic distri- 
bution of eccentricity would result, in which eccentricity 
would be small near the star, possibly large in the win- 
dow where 7 e > 0, and again small in the exterior region 
where j e ~ 0. The P or b range for large eccentricity may 
differ from one system to the next, as it depends on fac- 
tors such as stellar flux and planet mass, radius and tidal 
Q' p , and heating depth . This effect may explain the 
non-monatonic distribution of observed eccentricities for 
hot Jupiters. 

7. SIMULTANEOUS SPIN, ORBIT AND THERMAL 
EQUILIBRIUM 

In this section we self-consistently solve for the values 
of -Pspin, R P and e which satisfy 
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(73) 

Rather than choosing fixed values of F cmc independent of 
latitude, we follow appendix [Bl to compute F corc (9) by re- 
quiring that the radiative atmosphere match smoothly to 
an adiabatic interior with the correct M p and R p . Given 
this background solution for the atmosphere, we solve the 
boundary value problem for the linearized perturbations 
(see §3.5p . This includes diffusion explicitly, and inter- 
polates between the high and low frequency limits. The 
perturbations are then used to compute the quadrupole 
moments by integration over column feq |34[) and latitude 
(eqj2H|. Enou gh terms are used in the frequency sums in 
cq lAl71 IA18I and IA21I to ensure convergence. Newton's 
method is used find equilibrium values given a guess at 
the solution. We choose a small starting value for eccen- 
tricity, e = 0.01. If 7 e > 0, the eccentricity grows until it 
reaches a finite equilibrium value. If -f e < at e = 0.01, 
we treat e = as a good equilibrium solution. 

Results are shown in figure [14] for M p = Mj, = 
10 -2 cm 2 g _1 and Q' p = 10 6 . First we discuss the super- 
synchronous solution in the left panel. Close to the star, 
the gravitational tide dominates, resulting in nearly syn- 
chronous spin and circular orbit. The tidal heating rate 
due to asynchronous spin grows rapidly toward the star, 
leading to an increase in radius. The numerical values 
can be understood by comparing figures [TTJ and 1121 In 
between P or b — 4 — 6 days, 7 e > and e = is not a 
stable equilibrium. We find finite equilibrium eccentric- 
ity with values e = 0.0 — 0.2. Finite eccentricity does 
not seem to be as important as asynchronous spin in the 
tidal heating rate since no feature appears between 4-6 
days for R p . The eccentricity equilibria in figure [T4l are 
stable. For the e = equilibria, de/dt < implies sta- 
bility. In the window of instability, de/dt > for e < e cq 
while de/dt < for e > e cq . 

The timescale to return to spin equilibrium for a small 
deviation from equilibrium is 
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Fig. 14. — Equilibrium values of P sp j n , R p and e for M p = Mj, ft* 



io- 



g" 1 and Q' = 10 6 . 



super- (sub-) synchronous solution. In each panel, the solid line is for the left axis showing e, R p and P SI 
right axis showing the timescalcs to attain equilibrium. 

where N oq = \Neq | = |iVe q | is the value of the op- 



The left (right) plot shows the 
, , and the dashed line is for the 



posing torques in equilibrium. This timcscalc is shorter 
than the age by a factor ~ 10 3 for P or b < 1 week. The 
core cooling time in ea !61l is comparable to or shorter 
than the Gyr age for P or b < 10 days. The timescale to 
reach eccentricity equilibrium is 

/dine, 



dt 

(GT) 



(75) 



eval- 



where d In e eq /dt — |dlnee q /dt\ 

uated at equilibrium. The eccentricity evolution time is 
also comparable to or shorter than the age inside 1 week. 
Close to the star, we expect all three quantities to have 
equilibrium values, while further away there may not be 
sufficient time for R p and e to attain equilibrium values. 

The sub-synchronous solution is shown in the right 
hand panel of figure fT4l While the radius is similar to the 
that of the super-synchronous solution, the spin and ec- 
centricity differ. The (m, k) — (2, 1) resonance is weaker 
than the (m, k) — (2, 3) resonance due to the prefactor 
in ea l691 and eccentricity is driven less strongly than the 
super-synchronous case. At the evenly spaced orbital pe- 
riods used, we did not get close enough to the resonance 
for the driving to be significant. The spin period be- 
comes large for long orbital period, and will eventually 
turn into retrograde spin. 

Although not plotted in figure fT4l we found additional 
finite eccentricity solutions for P or b > 6 days, where 
e = is also a solution. We have not thoroughly in- 
vestigated large e equilibria which cannot be found from 
a small starting guess for e. We leave this for a future 
investigation, but note that many equilibria may be pos- 
sible due to higher order thermal tide resonances that 
become important for large e. The ability to trap the 
orbit into one of these resonances during evolution will 
require a more detailed numerical investigation. 

We now discuss how the results in figure [14] would 
change for different parameters. From figure [TT| increas- 



into balance requires that higher mass planets must have 
smaller radii, and vice versa for lower mass planets. For 
smaller Q' (larger gravitational tide dissipation rate), 
the planet must be further from the star before signifi- 
cant asynchronous spin develops. This pushes the strong 
resonance driving eccentricity further from the star (see 
figure fT3|). with a drop in the magnitude of the driving 
rate. Increasing F*(a) and e, and decreasing k*, acts to 
increase the strength of the thermal tide, pushing the 
spin away from synchronous and the region of 7 e > 
closer to the star where 7 e is larger. 

Our calculations are successful from the point of view 
that we have identified a powerful internal heat source 
which may explain the large radii observed for the tran- 
siting planets. However, our theory predicts a correlation 
of radius with orbital period, and nonzero eccentricity 
over a small range in orbital period. The plots of R p and 
e versus P or b for the observed transiting planets do not 
clearly show the simple behavior found in figure [Ml We 
speculate that part of the variance is due to differences 
in planet mass, stellar flux and heating depth («;*) as- 
sociated with the TiO/VO transition. Furthermore, our 
radii are upper limits, as we do not include a core. We 
intend to make more detailed comparison to data in a 
future investigation. 

8. DIFFERENTIAL ROTATION 

In this section we discuss the wind profiles induced by 
the thermal tide torque and estimate dissipation rates 
due to differential rotation. 

The time-averaged force per unit volume on a fluid el- 
ement due to the stellar tidal force acting on the thermal 
tide is 

f(x) = -(6pVU) t =p&VU)t. 



M leads to 



mg ^ p 
bring heating 
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larger cooling rates 
Rt) and cooling (oc 



at fixed R v 
M p R^) 



To 
back 



T 

In the high frequency limit, for heating function 
exp(— K-kU) we find a zonal component 



'. pn 2 R p 



(CI - n)C p T 



(76) 
Se = 

(77) 
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The ratio of the zonal acceleration = /^/pon a fluid 
element to the planet's gravity is then 



n 2 R Ti 



■10' 



(at 



th) 



M*_ (Rp 



(crtthY 



( M* \ ( 100 R p 



^) (atth)- 1 . (78) 



Vio 3 A/ p y v 

In the absence of friction and Coriolis forces, the torqued 
layers would accelerate to the sound speed on a timescale 
of days. Coriolis forces will act to bend the zonal winds 
to produce circulation patterns with characteristic veloc- 
ities w c ; r ~ a^/fl cos 9 sufficiently far from the equator. 
More importantly, the thermal tide force deposits net 
angular momentum. In the likely event that opposing 
gravitational tidal torques are exerted in a deeper layer, 
angular momentum transport is required to connect the 
source, due to the thermal tide, and the sink due to the 
gravitational tide. 

As a simplest possible model, we ignore the Coriolis 
force and consider a turbu lent boundary layer model (e.g. 
lLandau fe; Lifshitz~lll987f) that is forced by the thermal 
and gravitational tidal torques. The thermal tide exerts 
a stress 

T = jdzU = Jdyn*R PJ ^^ 

on the tidally forced layer. We treat this stress as applied 
in a thin layer at z = Z\, and the gravitational tide exerts 
an opposite stress at z = z-i < Z\, deeper in the layer. 
The horizontal momentum equation for our toy model 
with horizontal velocity v and kinematic viscosity v is 

d 2 v 

= P v—+t5(z-z 1 )-t5(z-z 2 ). (80) 

Setting the shear to be zero outsize zi < z < z±, we find 
the shear in between the two forced layers to be 

dz pv 

Turbulent velocity fluctuations of size 

1/2 



>(T/p)^=40 
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~ 5 g cm 3 



1/2 



1/2 



4 days 
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(82) 



Q-nJ \10 10 cm J \ P orb 

transport momentum. We take a characteristic vertical 
distance of a scale height H, so that the viscosity v ~ 
Hv+. Eq[HI]then becomes dv/dz ~ v+/H. The shearing 
implies an energy dissipation rate 

2 



-Eshcar ~ ^Rt 



dzpv ( ^ 
dz 



= 8 x 10 27 erg s" 
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(83) 



This estimate for dissipation due to shearing is compa- 
rable to the dissipation rate of the gravitational tide in 
eqU 



9. SUMMARY AND CONCLUSIONS 

We investigated the role of thermal tides in short- 
orbital period gas giant exoplanets, the hot Jupiters. Our 
investigation consisted of two parts. We computed the 
time-dependent temperature and density perturbations 
in radiative atmospheres. We then applied these results 
for the thermal tide, and the Darwin-Hut theory of the 
gravitational tide, to compute spin, radius and eccentric- 
ity evolution. 

Pre vious investigations (e.g. Ilngersoll &; Dobrovolskisl 
[19781) for the case of Venus considered a shallow at- 
mosphere bounded by ground at the base. In consid- 
ering the thermal tide in hot Jupiters, we used radia- 
tive heat transport, and allowed incoming stellar radia- 
tion to be fully absorbed in the atmosphere. We com- 
puted the temperature perturbations in the atmosphere 
resulting from time-dependent insolation in order to esti- 
mate the induced thermal tide quadrupole moment. As 
shown in in the context of Venus, in the high frequency 
limit, the thermal tidal bulge leads the subsolar point 
(|Ingersoll fc Dobrovolskisl[l978D . In considering deep at- 
mospheres, we include the effects of thermal diffusion, 
and found that the thermal tidal bulge still leads the 
substellar point for low tidal forcing frequencies, when 
heat can diffuse deep into the atmosphere. This result 
for the phase lead (for density, lag for temperature) is 
simplest to see in the results of nonlinear simulations in 
figures [2] - [5] Downward going diffusion waves at the 
base of the heated layer are discernible in figures [2] and 
[5] for the low frequency limit. The effective column and 
quadrupole moment are summarized in analytic formula 
in the high and low frequency limits in eqEJ |43l [49] and 

G2B 

In appendix [A[ we derive the evolutionary equations 
for the planet's spin and orbit, as well as tidal heat- 
ing rates for arbitrary quadrupole moments. These for- 
mula as used with the thermal tide quadrupole to com- 
pute equilibrium spin rates (SJ), equilibrium tidal heat- 
ing rates the growth/damping rate of eccentricity 
(Q, and simultaneous spin, orbital and thermal equilib- 
rium (<J7]). 

We find, contrary to the widely held assumption, that 
for reasonable heating depths and tidal Q, hot Jupiters 
are far from a synchronous spin state (see figure I10p . 
At long orbital periods, the thermal tide torque dom- 
inates and the planet is highly asynchronous, while at 
smaller radii the gravitational tidal torque is larger, 
and the planet approaches a synchronous state. Direct 
measurements of the planetary rotation rate could con- 
firm asynchronous rotation, although this appears dif- 
ficult in the near future. Possible methods include de- 
tection of centrifugal flattenin g using high photom etric 
precision transit observations (jSeager fc Huil 120021) . and 
Doppler shifted ab sorption lines for transiting planets 
([Spiegel et al.ll2007h . 

A more promising route to observationally constrain 
rotation rates at the moment may be through the indi- 
rect effect of planetary rotation on zonal winds, and the 
resulting change in phase of maximum thermal emission 
in Spitzcr infrared observations. We have found, for uni- 
form rotation, that flux perturbations are smaller in mag- 
nitude and have longer lag times in the high frequency 
limit than in the low frequency limit (see figure [4|) . The 
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strength of zonal winds is e xpected to decrease for more 
rapid rotation rate as well (jShowman et~aT1l2008h . Hy- 
drodynamic simulations of global circulation already con- 
tain the thermal tide quadrupole moment. "Turning on" 
the stellar tidal force in these simulations will torque the 
surface layer as described in $8j leading to spinup of the 
atmosphere, as well as momentum-driven circulation pat- 
tens in additio n to the usual the rmally-driven circulation. 
As argued by iGoodmanl ((2008), a source of dissipation 
such as that due to the gravitational tide or viscous fric- 
tion is necessary to damp kinetic energy added to the 
atmosphere, and may affect steady-state wind speeds. 

In Sj5j gravitational tidal dissipation due to asyn- 
chronous spin is shown to yield large steady state heat- 
ing rates, sufficient in magnitude to power the large ob- 
served radii over the range R p = (1.2 — 1.8) x Rj if the 
heat is deposited in the convective core. $5] shows that 
vertical shearing in the wind profiles between the layers 
torqued by the thermal and gravitational tides may also 
give heating rates large enough to alter the atmospheric 
temper ature profiles, and perhap s slow the cooling of the 
planet ( Bodenheimer et al.ll200lf ). 

In Sj6] we show that the thermal tide is capable of 
driving eccentricity. We compute the growth rate of 
eccentricity, in the limit of small eccentricity, finding 
narrow windows in orbital period in which the eccen- 
tricity is driven, while outside these windows it is even 
more strongly damped than for the gravitational tide 
alone. There are t wo co nsequences. As pointed out by 
iMatsumura et all (|2008l ). a puzzle exists as to why all 
hot Jupiter orbits are not highly circular. It is observed 
that ~ 25% of hot Jupi ters have eccentricities e > 0.1. 
IMatsumura et all (|2008h argue that this must imply the 
planets with eccentric orbits have far larger Q than plan- 
ets with circular orbits. Including the thermal tide, we 
find that it is expected that narrow windows in orbital 



period should exist in which eccentricity is driven to large 
values. This may explain the puzzling mix of zero and 
nonzero eccentricities in the observed planets. As a corol- 
lary, one cannot simply constrain the gravitational tidal 
Q, since it's effects are mixed in with the thermal tide 
effect. In some ranges of orbital period they combine to 
cause more rapid circularization, while in other orbital 
period ranges they partially cancel. A more detailed 
study including both effects is required to place mean- 
ingful constraints on Q. 

S}7| shows that the timescales to achieve spin, orbital 
and thermal equilibrium are small or comparable to the 
Gyr ages of observed transiting planets. Simultaneous 
equilibria for spin period, radius and eccentricity are 
found. Radii increase strongly toward the star, mainly 
due to asynchronous spin, inside orbital periods of 1-2 
weeks. Asynchronous spin again increases to longer or- 
bital periods. Nonzero equilibrium values of eccentricity 
are found in the region where e = is unstable. 

In summary, this initial investigation has found that 
thermal tides have a significant impact on hot Jupiter 
rotation rates, eccentricities and thermal state, and that 
seemingly unrelated observation puzzles, such as large 
radii and nonzero eccentricities, may have a natural ex- 
planation within this model. 
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APPENDIX 

SPIN, ORBIT AND THERMAL EVOLUTION RATES FOR ECCENTRIC ORBITS 

The interaction Hamiltonian coupling the density field of the planet to the gravitational tidal potential of the star 
is given by 



H = J d 3 xp(x,t)U(x,t) 



(Al) 



where the integral is over the perturbed body of the planet, p(x,t) is the time-dependent density in the planet and 
U(x, t) the tidal potential. To evaluate this expression, we expand the tidal potential in spherical harmonics 

U(x,t) = -GM i <J2^r L T (tjzt) r^frA Wm(M) (A2) 



Cm 



2£+l Ifl^ 1 



and define the time-dependent multipole moments of the planet 

Qirn(t)= I d 3 X r%* m (e^)p(x,t). 



(A3) 

Since r i p{x,t) is a real quantity, and Yg m = (— l) m Y^_ m , the mome nts m ust satisfy Q* lm = Qt,- m {— l) m . By defining 
Wim = [4%/(2£+ l)]Ye m (n/2, 0), the interaction Hamiltonian in eg I All may be conveniently expressed in terms of a 
sum over the tidal (spherical) harmonics 

ff(A*) = -GM*5^ m G2 m (i)_2_ (A4) 



The radial and tangential accelerations on the relative motion about the center of mass are then 

1 OH 
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wher e ij, = MyM+j (M v + M+) ~ M p is the reduced mass. The change in orbital energy and angular momentum are 
then (jMurrav fc DermottlH999l ) 

£ orb =/i (Dan +D$a$) (A7) 
i or b = ^Da®. (A8) 
Since i? rb = —GM p M±/2a and L or b = G{M P + M+)a(l — e 2 ), the changes in semi-major axis and eccentricity are 

a 2 



a navT" 



••2 



ee a Da$ 



[e sin $a_o + (1 + e cos $)a$] (A9) 



1-e 2 2a na 2 ^! - 

The torque on the planet is 



dU(x,t) f 3 ^dU(x,t) 



(A10) 



W = -y d d xp{x,t) K ' ' = J d 6 xp(x,t) ^ = -L orb (All) 

since the longitudes appear exclusively in the combination — $. This explicitly shows angular momentum is conserved 
over the orbit plus planet. 

Given the quadrupole moments Qe m (t), eq lA91 IA10I and lAlll can be integrated to find the changes in a, e and Q. 
To isolate the secular evolution, we treat the orbit as Keplerian with elements a and e that vary on timescales much 
longer than the orbital period. We then expand all quantities in a Fourier series in time. Let 

£+1 00 

M e -*n* = X l k m {e)e~ lknt (A12) 

k— — OG 

where the Hansen coefficients X^ n {e) are defined as (|Murrav fc Dermottl [l999) 



~5 mk + -[{1 + 1 + 2m)6 k ,m+i + (I + 1 - 2m)Sk, m -i] + 0(e 2 ). (A13) 
The Hansen coefficients are real, and satisfy Xj^ m = X^" 1 . The multipole moments are similarly expanded as 

oo 

Qim{t)= Qt™ke- iknt (A14) 

fc— — oo 

where 

/-2sr/ra 
n I • Jknt, 



Qe mk = 7T dt e QUt), (A15) 

27T Jo 

which satisfy Q L - m .-k = (-l) m Q; m 

Plugging eq |A"l2l and [Alii into eq. \M1 lAlOl and EH using the reality conditions relating (m, k) to (— m, — k), and 
and isolating the secular terms, we find 



a E, 



orb 



= n £ f (-2*) (A16) 



Imk 

,2, 



= n£ Wm xr (,»-*vT^) (A17) 



A ,mV £ W, m Xf" ( ^ftp ) M . (A18) 

£mfe ^ pfl 



The change in spin energy of the planet is 

E spin = QN = /mVfl Yl WtmXi™ ( ^ff^ ) (m) . (A19) 



The energy change in the orbit and spin is given by 



E mh + E spin = -/m 2 a 2 £ ( Im ^J fc) ) - mfi) . (A20) 



fnik 
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Here (Jkm = kn — m£l is the forcing frequency in the frame corotating with the planet. 

In the case of gravitational tides whose phase lag is due to dissipation, energy is taken out of the spin or orbit and 
deposited as heat in the planet, so E or b + E sp [ n < 0. In this case we define the heating rate of the planet as 



Keat - - (^orb + ^spinj = ^ 0? ^ W tm X^ | g^f ' j (kn - mQ) (A21) 



('ink 



This heating rate will be used for the thermal evolution of the planet. In the case of the thermal tide, work was needed 
to move the fluid against the stellar tidal force, so that E OT b + E spin > 0. Ultimately this energy is derived from the 
stellar radiation field. In this case we define the rate of work done on the atmosphere due to heating from insolation 
as 

E ( Z2 = (^orb + E spi n) ^ = -**nV £ W gm X l k m fe^LI) ) ( kn _ TO ft) . (A22) 

Imk \ p(1 J 

Since we do not include zonal motions in the equations solved in this paper, we cannot explicitly compute this work 
done against gravity. Our assumption is that the constant pressure approximation adequately represents the physics 
of heating by insolation. We will include fluid motions and compute this work in a future investigation. 

In equilibrium, the work done on the atmosphere first goes into spin or orbital energy, and then is dissipated as heat 
in the planet. Before equilibrium is reached, energy can be either lost or gained from the spin and orbit. 

The secular equations depend on the component of the quadrupole moments out of phase with the tidal acceleration, 
hence the imaginary components are required in the cq |A16[ IA17I and IA18I While we perform the sums over £, m and 
k until convergence is achieved for our numerical work, it is useful to have analytic limits to compare to. The torque 
on the planet for a circular orbit is dominated by the semi-diurnal term with \m\ = 2 and forcing frequency 2(n — fi). 
The torque at this order is 

37rV /2 fM p + AU 



N = A [To ) { Ml ) n Im ( Q222 ) ■ ( A23 ) 

The sign in eq |A23l is such that density perturbations leading the heating tend to torque the planet away from 
synchronous rotation, and vice versa. A schematic drawing for the relevant semi-diurnal tide (m = 2) is shown in fig. 

m 

At lowest order in e, the circularization rate is dominated by the m — harmonic with frequency n, as well as the 
|m| = 2 harmonics with frequencies 2(n — f2), 3n — 20, and n — 2Q. The circularization rate to lowest order in e is 

P:= \/f aT^ 2 " ^ eIm ( 2 222) + 3 Im(Q 201 )-7y|lm(Q 223 )-y|lm(Q 221 ) . (A24) 

Eccentricity may be either driven or damped depending on the relative sizes and signs of each term. 

MASS-RADIUS RELATION AND CORE COOLING LUMINOSITY 

Our model consists of an adiabatic interior connected to the "background" atmosphere structure described in i)3.5l 
We do not include a solid core at the center, which would decrease the radius. 
For the adiabatic interior we solve 

— =Anr 2 p (Bl) 
dr 

dr r z 

S(P, T) = S^c = constant. (B3) 

We use the equation of state from lSaumon et all (|1995l ) with 70% hydrogen and 30% helium by mass. Given the core 
entropy Score an d central pressure P c , the equations are integrated out to a reference pressure P rc f = 10 4 dyne cm -2 , 
yielding a mass M p , radius R p , and surface gravity g — GM p /R 2 . Note that the computed radii are at the fiducial 
pressure, not the self-consistent radius at the optical or infrared photosphere including transit geometry effects. The 
error in the radius will be of order a few scale heights, which amounts to of order a percent of the planetary radius. 

Given the solution for the adiabatic core, we integrate eq[29] and [30] inward from the top of the atmosphere. In 
addition to "given" parameters F»(a), e, 9 and «;*, a guess for F CQte must be chosen. We integrate inward until the 
Schwarszchild criterion is violated, yielding the entropy of the radiative atmosphere, Satm, at the radiative-convective 
boundary. In order that the radiative atmosphere matches smoothly onto the adiabatic interior, we adjust -Fcorc until 
Satm = Score- This yields the "background" radiative atmosphere, T versus P, used to find the perturbations in !j3.51 
The cooling luminosity of the core is found by integrating F cmc over latitude 

PIT 

icore = 2vri? 2 / dBwa.9 F coie {6). (B4) 
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An example of L cmc is given in figure 111! 

This model for the mass-radius relation and cooling luminosity is a modified version lArras & Bildstenl (|2006| ). Then- 
upper boundary condition assumed the existence of an isothermal region just below where the stellar radiation is 
absorbed. They assigned this regio n the temperature Td ee p, and treated this as a parameter of the model. Here we 
improve on lArras fc Bildstenl (|2006f ) by calculating, rather than parametrizing, the temperature at the atmosphere. 
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